Barrier crossing to the small Holstein polaron regime 
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We investigate the dimensionality effects of the Holstein polaron from the fully quantum regime, where the 
crossover between large and small polaron solutions is known to be continuous in all dimensions, into the 
limit described by the semiclassical Discrete Nonlinear Schrodinger (DNLS) Equation, where the crossover 
is continuous in ID but discontinuous in higher dimensions. We use exact numerics on one hand and a two 
variable parametrization of the Toyozawa ansatz on the other in order to probe the crossover region in all 
parameter regimes. We find that a barrier appears also in ID separating the two types of solutions, seemingly 
in contradiction to the common paradigm for the DNLS according to which the crossover is barrier-free. We 
quantify the polaron behavior in the crossover region as a function of the exciton overlap and find that the barrier 
remains small in ID and tunnelling through it is not rate-limiting. 



The Holstein model has been used widely to describe trans- 
port properties of electrons or excitons in diverse systems, 
ranging from molecular crystals ^ to strongly correlated 
electron-phonon systems 131, interface charge localization in 
alkane layers iH, organic transistors jstl, proteins and 
protein-like crystals S S Hi as well as DNA lH]. In the 
context of the model, the electronic degrees of freedom are 
coupled to dispersionless phonons leading to the formation of 
large, i.e. much larger than the lattice spacing, or small, viz. 
essentially single-site, polarons. The relative magnitudes of 
three parameters of the theory, viz. electronic overlap inte- 
gral, phonon energy and exciton-phonon coupling determine 
the specific polaron properties. Outstanding issues have been 
the nature of polarons in different regimes, their transitions as 
well as the precise onset of selftrapping. The Discrete Non- 
linear Schrodinger (DNLS) equation, that is a semiclassical 
approximation to the Holstein model, gives a continuous tran- 
sition in ID and a discontinuous transition in higher dimen- 
sions ifiol [Till , while most variational calculations, in partic- 
ular those based on the Toyozawa wavefunction, show a dis- 
continuity in all dimensions |ll[ll[ll[ll[ll[Il[Il- On 
the other hand, it has been shown recently with the help of a 
numerically exact solution of the full-quantum Holstein po- 
laron that the cross-over is continuous in ID lfl9ll as well as 
in higher dimensions ll20l I21I1 . The non-existence of phase 
transitions in polaron systems can also be proven on very gen- 
eral grounds 1I22I1 . In the present paper, we address the nature 
of the polaron transitions using exact numerics for the fully 
quantum Holstein model on one hand, as well as a simple, 
physically motivated two parameter variational minimization 
on the other hand. The comparison of the two approaches 
resolves the issue of the large to small polaron transition, pro- 
duces quantitative information on the nature of the selftrap- 
ping transition and finalizes pending issues regarding deficien- 
cies of standard approximations. 

The Holstein Hamiltonian reads: 

H = Hex+Hpl,+Hex,ph 

He.. = -JJ^[b/Bj+,+b;Bj^,) 
j 
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He...pi, = -xLB/Bj{b/ + bj) (1) 
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where b^j (bj) and fij (Bj) destroy (create) a phonon and an 
electron (or exciton) at site j, respectively, and /, h(0, % are the 
excitonic overlap, phonon energy and exciton-phonon cou- 
pling, respectively. In what follows we use two dimension- 
less parameters, i.e. the exciton coupling J/h(0 and exciton- 
phonon coupling x/h(£) while all energies are given in units 
of one phonon-quantum fico. We also restrict our study to the 
wavenumber k ^0 case since this will be the polaron ground 
state, and because it is most relevant in optical spectroscopy 
(k = Q selection rule). 

We first present numerically exact results in all three dimen- 
sions using the approach of Trugman and coworkers lfl9ll20ll . 
Figure [1] top row, shows the reduced phonon density at the 
site of the exciton as a function of the phonon-coordinate 
^0 = (^0 + ^0), with phonon coordinates at sites different from 
the exciton site traced out. We note that the results are qualita- 
tively similar in all dimensions: The reduced phonon density 
is centered at q'o ~ X for small exciton couplings J (i.e. a small 
polaron), similar to the 7 = case, for which an exact (analyt- 
ical) solution exists revealing with ^0 = X tSI- For sufficiently 
large exciton couplings J, the phonon displacement shifts to 
^0 ~ 0, i.e. essentially a free exciton without any phonon- 
displacement. We observe that the transition between the two 
regimes is smoother in ID than it is in 2D and 3D, however, 
the essential point is that free exciton and small polaron solu- 
tions coexist in a certain parameter range, also in ID. The two 
states gradually change their relative weights, but hardly their 
character, as the exciton coupling J is increased. If we inter- 
pret the reduced phonon density as one that originates from 
an effective potential V{qo), then the coexistence of two solu- 
tions hints to the presence of a barrier separating them. This 
conclusion is in disagreement with the semiclassical DNLS 
case ifiol [Till , which does not reveal any barrier in ID but 
rather a gradual transition. 

The phase diagrams in Fig.[T] bottom row, summarize these 
results where we plot the expectation value of the mean dis- 
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FIG. 1: (Color online) Results from a numerically exact diagonalization of the Holstein Hamiltonian in ID (left), 2D (middle) and 3D (right). 
Top row (a-c): Phonon density of the polaron ground state as a function of the phonon-coordinate q^). Shown is the reduced phonon density at 
the site of the exciton. The exciton-phonon coupling was set to x = 4.0. Bottom row (d-f): Phase diagram, plotting the expectation value of 
the phonon coordinate (^o) in units of % for the polaron ground state. Contour lines are equidistant with spacing 0.1. In red (light gray) are the 
phase separation lines from Eq.[2] whereas in blue (dark gray) are those obtained from equalling the two variational minima of the Toyozawa 
wavefunction (see Fig.|2j. 



placement of the phonon coordinate {qo) in units of % (which 
is the maximum value obtained in the 7 = case). For small 
couplings X and /, the result is qualitatively very similar in all 
dimensions with a continuous transition between small and 
large polaron solution. Nevertheless, it is quite evident that 
the transition becomes increasingly more abrupt in 2D and 
3D for large couplings (i.e. as we approach the regime of the 
DNLS), whereas the cross-over stays smooth in the ID case. 
The phase separation lines can be fitted extremely well to a 
generic relationship (Fig.[T] bottom row, red lines): 



aJ 



(2) 



with a = 1 in ID, a = 3.34 in 2D and a = 5.41 in 3D, re- 
spectively. The ID value stems from Lindenberg et al. em- 
pirical relationship [23J, whereas those for 2D and 3D co- 
incide with the critical couplings in the DNLS limit (when 
\/aJ ^ 1 /a) ifTHl . Hence, Eq. (O does give the correct semi- 
classical limits expected from the DNLS analysis. 

Although exact, the results of Fig. [T] originate from a di- 
agonalization of a huge matrix and provide relatively lit- 
tle physical insight. We therefore employ an approximate, 
yet more intuitive trial function, the Toyozawa wavefunc- 
tion lfl2i [13I [Til [itIi , and use the numerically exact solution 
as an accuracy check. The Toyozawa trial function is written 
as a Bloch state: 



1^) = 



iki 
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with 



(4) 



where all |\|//) = |\|/) are identical, lO)^ is an exciton vacuum 
state and \q)j a coherent phonon state: 



(5) 



The Toyozawa wavefunction function can be considered a 
Bloch-extension of symmetry breaking soliton solutions 13, 
where both excitons and phonons are dressing one common 
trapping site /. In ID, the norm of the wavefunction \\\t) 
isllJ: 



(6) 



and the expectation value of the Hamiltonian: 
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FIG. 2: Variational energy of the Toyozawa wavefunction as a func- 
tion of qo and in ID (left), 2D (middle) and 3D (right). Exci- 
ton coupling is J = 9/D and the exciton-phonon coupling x = X™' 
has been calculated so that both variational minima are equally deep. 
Contour line spacings are 0.05 in (a) and 0.1 in (b) and (c). 



In more than ID (dimension D), the indices / and j are re- 
placed by and {jxJyiJz)}' respectively, in either 
sum/product. Furthermore, the exciton term is repeated D 
times with the one-site shift (i.e. the ^/_,+i-term) in either 
direction (for isotropic exciton coupling, as considered here, 
this can be reduced to one exciton term with prefactor 2DJ). 

Zhao et al. ifTill have minimized the energy of the Toyozawa 
wavefunction by effectively varying all parameters {«,} and 
{qi} (the calculation was done in momentum space), still leav- 
ing us with a multi-parameter solution. In contrast, we make 
an exponential ansatz, in analogy to Refs. fill [Tsl [itIi : 

I 'I 

qi = qoy\% (8) 

that contains only three parameters. Minimizing Eq. [T] re- 
veals qi as one condition ifisi [Till which allows us to 
eliminate one of these parameters leading to r|p/, = ( ^/^ — 
\/qo)/{\/X+ v^)- Hence, we can express the variational 
energy as a function of effectively two parameters, and rj^. 

Figure |2] shows the variational energy of the trial function 
for 7 = 9 /D. We observe the emergence of two variational 
minima in all dimensions with a barrier separating them. The 
barrier is very shallow in ID (« 0.1) while it is much more 
pronounced in 2D or 3D. This is the barrier which is respon- 
sible for the double-peak structure of the wavefunction in the 
numerically exact solution, which also exists in ID (Fig.[T^). 
Venzl and Fischer ifisll have shown a very similar variational 
energy surface in ID (as a function of r|« and r|p/,, eliminat- 
ing ^o). however, did not address the double-minimum due to 
a choice of too small value of te exciton coupling, viz. 7 = 2. 

For the chosen values of the exciton-phonon coupling the 
two minima are equally deep. Using the degeneracy of 
the variational minima as a criterion for the critical exciton- 
phonon coupling x'^"'' that induces the large-to-small polaron 
transition, we find that the numerically obtained line Xc"' i-^) 
follows very closely the relationship of Eq. |2] (compare blue 
and red lines in Fig.[TJl-f). This near identity of x"" and Eq.|2] 
starts to deviate in ID for very large exciton overlaps J > 100, 
indicating that Eq. |2] is a lower order expansion of a more 
complex functional relationship. 

Figure [3] shows the energy of the lowest energy eigenstate, 
as deduced from the numerically exact solution (black line). 
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FIG. 3: (Color online) Energies of the lowest energy state for 7 = 9 
and the exciton-phonon coupling % varied through the cross-over re- 
gion in ID. Black line: numerically exact result. Red (light gray) 
lines: The two minima of the variational energy. Blue (dark gray) 
line: Energy corrected by tunnelling between the two components. 
The inset focuses into the crossing region and defines the terms in- 
troduced in Eq.lIOl 



and compares it with the two variational minima (where they 
exist, red lines). Outside the crossing region, the variational 
result agrees reasonably well with the exact energy - hence 
the two-parameter trial function Eq. [8] is quite good in this 
regime. However, both deviate significantly at the crossing 
point, indicating that the Toyozawa wavefunction per se does 
not capture the physics of the crossover correctly 11121 [Til [Tsll . 
In particular, the variational minimum, as well as all other 
properties of the trial wavefunction, would not be analytical 
at the cross-over point. This is simply due to the fact that the 
two variational minima are separated by a baiTier (Fig. |2^), 
and the solution of an overall minimization switches abruptly 
when one of these local minima decreases below the other 

Barisic has extended the variational space by considering 
a linear combination of two (or more) Toyozawa wavefunc- 
tions HI: 

|v) = q|V")+c..|V^'), (9) 

showing that this removes the discontinuity in ID. Fig. |2^ 
makes very clear why that is so: The two Toyozawa wave- 
functions will sit in either of two minima of the variational 
energy, and tunneling coupling between them will (a) lower 
the energy and thereby correct for approximately half of the 
mismatch to the numerically exact result, and (b) render the 
cross-over continuous (see Fig. [3] blue line) ll24ll 

In order to understand the nature of the barrier separating 
the small and large polaron solutions we address the question 
of how the crossover domain behaves in the limit J ^°o. Mix- 
ing of the large and the small polaron states will occur when- 
ever the tunnel coupling is on the order of the energy separa- 
tion between them, Ztimnd ^ —£'■'■'1. Hence, the quotient 
of the tunnel coupling divided by the difference of the deriva- 
tives of the two variational energies with respect to x defines 
a measure of the extension of the crossover region (Fig.|3] 
inset, shows pictorially ttunnei and A^): 
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FIG. 4: (Color online) (a) Tunnel coupling Etumiel and (b) extension 
of the cross-over region A/ as a function of exciton coupling J (with 
X = X™' ) in ID (red), 2D (blue) and 3D (green). 

Figure|4^, red line, shows that the tunnel coupling ttwmei stays 
essentially constant in ID for all 7's, however, the denomina- 
tor in Eq.[TO]decreases with J so that A% increases roughly lin- 
early for large enough / (Fig.|4]3, red line). Consequently, the 
cross-over remains continuous even for J ^ °o, leading to the 
correct ID result in the adiabatic limit. In contrast, in higher 
dimensions, both the tunnel coupling ttimnei as well as de- 
cay exponentially with / (Fig.|4] blue and green lines). Hence, 
the superposition of two Toyozawa wavefunctions, although 
predicting in principle a continuous crossover, leads (in 2D, 
3D) to an increasingly more abrupt transition for 7 ~> oo, just 
like the numerically exact result does (Fig.[T^,f). 

In conclusion, we studied the dependence of quantum po- 
laron regimes on lattice dimensionality employing exact nu- 
merics and a physically motivated variational ansatz based on 
the Toyozawa wavefunction. The exact numerical solution 
shows that the large-to-small polaron crossover is continuous 



in all dimensions although much more abrupt in two and three 
dimensions, especially as the parameters approach the adia- 
batic regime (Fig.[T^,f). In ID, the transition is smooth and 
remains continuous also for large /'s (Fig.[T}l). 

More physical insight into the polaron features is obtained 
from an approximate trial wavefunction, viz. the Toyozawa 
wavefunction of Eq. |4] that leads to variational results whose 
regime of validity has been analyzed extensively lfl4l [Tsll . 
When the latter is augmented with an exponential ansatz of 
Eq. |8] containing effectively only two variational parameters, 
it leads to a very intuitive and clear picture. Despite the sim- 
plicity of this approach, it provides results that coincide with 
the exact ones outside the crossover region, while, neverthe- 
less, predicting non-analytic properties at the transition point 
in all dimensions. The appearence of this discontinuity was 
pointed out previously IUtI I20I1 but it has been found in this 
work to be due to the presence of two, almost degenerate solu- 
tions coexisting in a certain parameter regime while separated 
through a barrier. Consequently, a superposition of two Toy- 
ozawa wavefunctions, as in Eq|9] reveals qualitatively correct 
results for the cross-over region in all dimensions and in all 
parameter regimes iHtIi . The resulting intuitive picture may 
be used to asses quantitatively the properties of the cross-over 
region as a function of the exciton coupling. One important 
outcome of this analysis is that a barrier separating large and 
small polaron solutions exists also in ID, in sharp contrast 
to the predictions of semiclassical theories. We nevertheless 
find that tunnelling through the barrier remains efficient in all 
parameter regimes (Fig. |4^, red line), since the barrier stays 
of the order of /ico and hence does not hamper motion in the 
DNLS limit where phonons are treated classically. The Hol- 
stein model is a fundamental "minimal" model for strongly 
interacting systems and thus our results may also have ram- 
ifications for more complex models such as for instance the 
Holstein-Hubbard model. 
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term stemming from |\|/'''), and the other, aj^i or q^-h from creasing the mixing matrix element on the cost of increasing 

|\|/("'). Furthermore, in the ID case, we have allowed the po- the diagonal matrix elements (energies), 

sitions of the two Toyozawa wavefunctions to deviate slightly 
from the exact minima during the minimization, thereby in- 



